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ABSTRACT 

Aims. Chemical networks used for models of interstellar clouds contain many reactions, some of them with poorly determined rate 
coefficients and/or products. In this work, we report a method for improving the predictions of molecular abundances using sensitivity 
methods and ab initio calculations. 

Metliods. Based on the chemical network osu.2003, we used two dilferent sensitivity methods to determine the most important 
reactions as a function of time for models of dense cold clouds. Of these reactions, we concentrated on those between C and Cj and 
between C and Cs, both for their effect on specific important species such as CO and for their general effect on large numbers of 
species. We then used ab initio and Idnetic methods to determine an improved rate coefficient for the former reaction and a new set of 
products, plus a slightly changed rate coefficient for the latter 

Results. Putting our new results in a pseudo-time-dependent model of cold dense clouds, we found that the abundances of many 
species are altered at early times, based on large changes in the abundances of CO and atomic C. We compared the effect of these new 
rate coefficients/products on the comparison with observed abundances and found that they shift the best agreement from 3 x 10"* yr 
to (1 - 3) X 10"^ yr 

Key words. Astrochemistry - Molecular processes - ISM: abundances-molecules-clouds 

1 . Introduction "Important" reactions can be defined in several ways. First, 

one can define important reactions for specific species; such 
0\ The chemistry in astrochemical models is driven by large net- reactions strongly affect the abundances of individual species 
^ works including thousands of reactions with poorly determined through changes in their rate coefficients for formation or de- 
rate coefficients. For several years now, astrochemists have been struction processes. One must be precise with this definition, 
working on the quantification of model uncertainties caused by because even though some reactions dominate the formation or 
uncertainties in rate coefficients. Using similar methods but dif- destruction of certain species, changes in t heir rate coe f ficient s 
ferent models and networks, Vasyunin et al. (2004, 2008) and within certain ranges can have little effect. lOuan et al.l (l2508h . 
■ - 'IWakelam etal.| ^005, 2006aD have studied the error propaga- for instance, showed that although molecular oxygen is formed 
tion of such uncertainties in various environments such as dense by the reaction O + OH ^ O2 + H, the rate coefficient of this 
clouds, diffuse gas, hot cores, and protoplanetary disks. This ap- reaction can be changed by several orders of magnitude without 
proach to quantify the model eiTor is useful quantitatively com- affecting the O2 abundance significandy. The reason is that the 
paring calculated and observed molecular abundances. In addi- reaction is also the main route of OH destruction. When the rate 
tion to this important aspect, these authors have developed sen- coefficient is decreased, the OH abundance increases so that the 
sitivity methods for identifying those reactions in complete net- rate of O9 formation stays constant. As a counterexample, the 
works that might be worth studying theoretically or experimen- reaction C + C3 ^ C4 + hv is very important to the abundance of 
tally because model results are very sensitive to their rate co- cO, as will be discussed below. Secondly, "important" reactions 
efficients and uncertainties. Identifying these reactions is criti- can be defined through the number of molecular abundances they 
cal because detailed experimental or theoretical studies of a sin- influence strongly. As an example, the cosmic-ray ionization of 
gle reaction, especially at low interstellar temperatures, can take molecular hydrogen is important for almost all of the species in 
years. Guiding chemists and physicist as to which reactions are the chemical models of dens e interstellar clouds (.Wakelam et all 
the most important for them to study is an efficient method for |2006bl : IVasvunin et alj|2008h . 

reducing error in chemical models. ^^-^ p^p^^^ p^^^^^^; ^ complete approach to improv- 

ing chemical models through identifying sensitive reactions fol- 

Send offprint requests to: V. Wakelam, e-mail: lowed by an increase in our knowledge of their rate coefficients 

wakelam@obs.u-bordeauxl.fr via theoretical means, and using the newly revised rate coeffi- 
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cients in model predictions. We direct our focus on two reac- 
tions: the radiative association reaction between C and C3 and 
the neutral-neutral reaction C -1- C5. We start, in Section 2, with 
a discussion of the methods used for our sensitivity analyses 
and the determination of the important reactions in cold dense 
clouds. The theoretical study of the C -n C3 and C -H C5 rate co- 
efficients is presented in Section 3 and the impact of new rate 
coefficients on the model predictions are studied in Section 4. 
Following our conclusion, an Appendix on other reactions of the 
C + Cn class is presented. 



2. Selection of important reactions 

In order to identify important reactions of both classes, we use 
two different sensitivity methods and compare the results of the 
two. Each sensitivity method is based on linear correlations but 
in the first one, only one reaction is modified at a time whereas in 
the second method all reactions are modified at the same time. 
The first one, which we will call "individual perturbations", is 
simply based on the study of the linear effect of individual rate 
coefficient perturbations on molec ular abundances. We have al - 
ready used this method previously (IWakelam et al.l2005ll2006al) . 
In the second method, we compute Pearson correlation coeffi- 
cients. 



2.1. Chemical model 

To represent cold dark clouds, we use a pseudo- time-dependent 
gas-phase model, which computes abundances as a function of 
time for a temperature of 10 K, a total hydrogen density of 
2x10'* cm"^ and a visual extinction of 10. The model is known as 
nahoon (http://www.obs.u-bordeauxl.fr/amor/VWakelam/). For 
comparison with previous studies, the elemental abundances 
are the O-rich low-metal elemental abundances (see table 1 in 
IWakelam & Herbstll2008l) . For this work, we used an older ver- 
sion of the osu database (osu.2003) since one goal of this work 
is to show how we corrected a mistake in this database. We focus 
on the specific time of 10^ yr, since this so-called "early time" 
normally results in the be st agreement betw een calculated abun- 
dances and observations dSmith et alj|2004l) . 



2.2. Sensitivity Methods 

The idea of the individual perturbation approach is to look at 
the individual effect of each reaction. With a network of 
reactions, we run the model times and modify only one rate 
coefficient each run. One then obtains the chemical evolution of 
the fractional abundance: -^y(f) for each species j and run / . Note 
that the index / can stand for the run or the reaction with varied 
rate coefficient. We then compare the abundance of each species 
^'jit) with the standard abundance Xj^it) obtained without 
changing the rate coefficients. In our previous uses of this ap- 
proac h to select importan t reactions for hot-core jWakelam et alj 
l2005h and dense-cloud ( IWakelam et alj l2006al) chemistry, we 
multiplied each rate coefficient by its uncertainty factor This 
approach assumes that the real uncertainty in the rate coefficient 
is not much larger than the considered uncertainty factor, which 
may not be true. Indeed, as we will show later, this approach 
can miss important reactions if the real rate coefficient error is 
much larger than the uncertainty factor For this reason, here we 
multiply each rate coefficient by the same factor of two. This 
factor is currently the default uncertainty used for most classes 
of reactions in current reaction networks ( http://www.udfa.net/, 



http : //www . physics . ohio- state . edu/~ eric/research. html). 
To quantify the effect of these variations in the reaction rate 
coefficients on the abundances, we compute a parameter labeled 
R'jit) for each species and each reaction, where 



(X'j(t)-Xf(t)) 
Xf(t) 



(1) 



As an example, a value of ^^(0 = 0.2 tells us that the multiplica- 
tion of the reaction / by a factor of two increases the abundance 
of species j by 20%. We also did these calculations by multiply- 
ing the rate coefficients by a factor of 10 or dividing them by 2 
but the list of important reactions is rather the same as for the 
case presented here. 

The second method is based on linear correlations but is less 
time-consuming. Here we run about 2000 different models I, ran- 
domly varying all rate coefficients within an uncertainty range 
of a factor of two, with both increases and decreases consid- 
ered. We then compute Pearson correlation coefficients for each 
species j and each reaction ; using the following relation: 



Y.'iX'jit)-Xjitm',-kd 
^il!(X'ft)-X]{t)fY!ik\-ki)^ 



p'jif) = 



(2) 



where Xj(t) and kj are the mean values of the abundance for 
species j and the rate coefficient /. Pearson correlation coeffi- 
cients indicate the strength and direction of a linear relationship 
between two variables by dividing the covariance of the two 
variables by the pr oducts of their standard deviations (see for 
instance iPress & Flannerv..l992) . 

2.3. Results of sensitivity calculations 

In this section, we show the results of both sensitivity methods 
to determine the most important reactions for the abundances of 
selected individual species and to determine the reactions that 
influence the abundances of the largest numbers of species. 

In Table [1] we show the most important reactions for the 
species CO, O2, C4H and HC3N at an integrated time of 10^ yr 
using the individual perturbation (R) and the correlation coeffi- 
cient (P) methods. In tables |2] and [3l we list the reactions that 
influence the most species using R and P, respectively, as a cri- 
terion. For each reaction, we count the number of species for 
which Roi P is either larger than 0.1 or smaller than -0.1. 

A glance at Table [T| shows that the important reactions gov- 
erning the abundances of individual species are generally sim- 
ilar for both methods. There are two main surprises, however. 
First, for the case of CO, the values of R are much smaller than 
those for P. Presumably, the explanation for this dichotomy is 
that most of the carbon goes into CO at that time. As a result, 
a change in the rate coefficients does not change the CO abun- 
dance, which is then controlled by the elemental abundance of 
C. This effect al so explains the small un certainty in the CO abun- 
dance found by IWakelam et al.l (l2006ah . Secondly, the most un- 
expected reactions, because they seem so unrelated to the listed 
species, are those between C and C3 and between C and C5. Note 
that the latter was assumed, erroneously, in osu.2003 to be a ra- 
diative association (see below). The importance of the C + C3 
association seems to be associated with the fact that C3 acts as a 
catalyst to convert C and O into CO via the association reaction 
between C and C3 followed by 



O + C4 — >C0 + C3, 



(3) 
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Table 1. Most important reactions for selected species at 10^ yr. 



Species 



CO 



O2 



C4H 



HC3N 



Reaction 



R 



C + C3 
C + C5 



> C4 + hv 
H+ + e- 

► Cs + hv 



0.009 
0.004 
0.005 



0.5 
0.3 
0.4 



H2 + ^ — 
C + C3 — 
H+ + CO 
C + C5 — 
C + O2 — 



H+ + e- 

> C4 + hv 

HCO+ + Ho 

> Cfi + hv 

> co + o 



1.2 
0.3 
-0.3 
0.2 
-0.2 



0.7 
0.2 
-0.2 
0.2 
-0.2 



O + C4H - 
He-Hf — 
C + C, — 
C + C3 — 
C4H^+ e- 



C3H + CO 
He+ + e- 
■Ce 
■ C4 

— ^ C4H + H 



-0.7 
0.3 
-0.2 
-0.2 
0.2 



-0.2 
0.2 
-0.07 
-0.08 
0.09 



He + ^ — 
H2 + f — 
H3O+ + e 
H3O+ + e 
H+ + C - 
C+ + H2 - 
C+ + HC3 
CH+ + H: 
C + C, - 



He+ + e- 

H+ +e- 



OH + H + H 
H2O + H 
CH+ + H2 
^ CH+ 
N — > C3H+ + CN 
CH+ 

C4 



-0.5 
0.5 
-0.2 
0.1 
0.1 
0.1 
-0.1 
0.1 
-0.1 



-0.5 
0.5 
-0.1 
0.2 
0.2 
0.2 
-0.1 
0.1 
-0.09 



Table 2. Reactions influencing the most species at 10^ yr using 
R 



Reaction 


% of species 


H2 + 4" — ^ H+ + e- 


53 


H+ + C — ^ CH+ + H2 


31 


He + 4- — > He+ + e" 


31 


C+ + H2 — ^ CH+ 


29 


H3O+ + e- — » H2O + H 


26 


HCO+ + C — ^ CH+ + CO 


24 


CH+ + H2 ^ CH+ 


22 


HCO+ + C3 — ^ C3H+ + CO 


16 


C + C3 ^ C4 


15 


H+ + CO ^ HCO+ + H2 


14 


C+ + e- — > C2 + C 


13 


H3O+ + e- OH + H + H 


13 


CH+ + CO ^ HCO+ + CH4 


13 



which regenerates the C3. A similar cycle happens with C5 if it 
undergoes radiative association with C. 

Considering the number of species affected, the most impor- 
tant reactions determined by both methods are quite similar al- 
though the order o f importance can differ. As already noticed 
by several authors (IXVakelam et al ]l2006aHVasvunin et alJl2004l: 



IWakelam et al.ll2006bl:IVasvunin et alj|2008l) . the H9 and He ion- 
ization rates by cosmic-rays appear to be important for many 
species since they are at the origin of formation and destruction 
pathways. Most of the other important reactions are ion-neutral 
normal and associative reactions and dissociative recombination 
reactions, such as those for H^O^ with electrons. Among the 



+ H2 ^ 

m + CO 



CH^, CHt + H2 



+ H,, Ht +0 



CH+, 



important reactions, 1 112 — ' ^^^2 
H+ -H C ^ CH+ + H2, H+ -H CO ^ HCO + 02, ..^ 
OH"^ + H2 have al so been identifi e d as im portant reactions for 
disk chemistry by IVasvunin et al.l (l2008l) . Most of these reac- 
tions are at the start of a chain of reactions that builds up hy- 
drocarbons. It should be emphasized that our sensitivity meth- 
ods need not tell us the details of the pathway of reaction cy- 
cles. If, for example, a destruction reaction leads to a loop of 



Table 3. Reactions influencing the most species at 10^ yr using 
P 



Reaction 


% of species 


Hj + ^ — ^ H+ + e- 


86 


He + f — > He+ + e" 


82 


H+ + — ^ 0H+ + H2 


54 


H3O+ + e- — ^ OH + H + H 


48 


H3 + CO — > HCO + H2 


43 


C + C5 Cs 


41 


C + C3 ^ C4 


41 


C+ -H H2 ^ CHJ 


41 


H3O+ + e- — » H2O + H 


37 


HCO+ + e- ^ CO + H 


34 


CH+ + H2 — ^ CHJ 


33 


HJ + C — ^ CH+ + H2 


32 


HCO+ + C — ^ CH+ + CO 


24 



subsequent reactions that reforms the molecule being destroyed 
with some non-zero efficiency, the sensitivity method might just 
decide that the initial destruction reaction is less important, at 
least for the species being destroyed. The problem of enhanced 
molecular abundances due to artificial loops has been minimized 
whenever possible in the construction of the network by includ- 
ing a number of different product channels. For example, if the 
major destruction route of a neutral species leads to a protonated 
ion, and the major destruction route of the ion is dissociative re- 
combination to form neutral molecules, the channel that loops 
back to the original neutral species plus a hydrogen atom will 
not be assumed to occur on more than 50% of the recombination 
reactions (unless, as is highly unusual, experimental results tell 
us the contrary). Let us look at the important reactions causing 
the destruction of HC3N in Table[T] One of them is the rather di- 
rect destruction via C^ ions. The products of this reaction do not 
lead back to HC3N in any direct fashion. If, on the other hand, 
there were a charge exchange reaction to form CsHN^, and this 
ion were hydrogenated by reaction with H2 to form, C3H2N^ 
+ H, the hydrogenated ion could react with electrons to form 
HC3N + H. In the network, however, the probability of this neu- 
tral product channel is set to approximately 50% so that the loop 
is not particularly efficient. 

The only neutral-neutral reactions that seem important for 
many species in Tables |2] and |3] are the (presumed) radiative as- 
sociative reactions C + C3 and C -i- C5. The rate coefficients for 
these rea ctions used i n osu. 2003 are 10""^ and 10"" cm-' s"' 
based on ISmith et aTl (l2004t) . without tempe rature dependence . 
These reactions have been looked at again bv lSmidi et al](l2006h 
and it was realized that the reaction between C and C5 has the 
normal exothermic products C3 -H C3. Given the importance of 
these reactions, we decided to study them theoretically in some 
detail. 



3. New calculations for C + C3 and C + C5 reactions 

3.1. C + C3 

At low collision energies, the reaction of C(^P) atoms with lin- 
ear C3 in its ground electronic state (X 'Sg ) leads only to adduct 
formation, C4, because the exit channel C2 + C2 is endothermic. 
The lowest energy pathway to produce C4 is the linear one, and 
we focus on this pathway here. Under interstellar low-density 
conditions, the only way to stabilize the C4 adduct is by emission 
of a photon. The formation of the linear C4 adduct is strongly 
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exothermic and spin-allowed to its triplet manifold of states, des- 
ignated ^€4: 

C(^P) +^ C3 -^^ C4 + hv. (4) 

There are a variety of values for this exothermicity. A high- 
temperature experiment yields a value of (5.21 +0.21) eV fr om 
the heats of formation of C4 and C3 (Gingeri ch et al.l ll994'). a 
PST calculation based on a photo dissociation expe riment yields 
a lower value of (4.84 +0.15) eV jChoi et alj|2000l) . while high- 
level ab initio theoretical values (including our own work) com- 
puted at the RCCSD and RCCSD(T) l evels range from 4.94 - 
4.98 eV without zero-point correction dMartin & Tavloil [19951) . 
The high-temperature experimental value has been used in the 
RRKM calculations discussed below. The rhombic isomer of C4 
in its ground singlet state has been calculated t o to b e nearly 
isoenergetic with the linear isomer (Masso et al] |2006l) . but the 
pathway from C and linear C3 to the first excited triplet state 
shows a high barrier at the CASSCF level of ab initio theory 
(see below). 

In the Coo.v symmetry, the three electronic states of C4 cor- 
relating to the three spin-orbit states of C(^P()_i_2) lead to two 
different potential energy surfaces corresponding to the funda- 



Table 4. Vdz/CASSCF optimized structural parameters 



mental state X 



Eg and the first excited 



'n^ state. (Note that 



the ground electronic state of C is of g symmetry.) The ener- 
getic and structural parameters of these states of C4 and the 
second excited ^n„ state are listed in Table H] These parame- 
ters and others such as the pathway to the rhombic isomer of C4 
have been calculated at the Restricted Hartree-Fock (RHF) and 
Complete Active Space Self Consistent F ield (CASSCF) level 
with the double-zeta cc-pVDZ basis set of Dunning (1989), us- 
ing MOLPRO 2000, a package of ab initio programs designed 



by H.-J. Werner, and P. J. Knowles (see http://www.molpro.net 



for more details). This level of theory was also used to deter- 
mine if the three electronic curves are attractive along the linear 
reaction pathway. In these calculations, which underestimate the 
exothermicity of reaction (|4|, the Coov (linear) approach of C to- 
ward the C3 axis, shown in Figure [1] yields no barrier for the 
ground C4(^Ep electronic state but shows a high barrier of 1.5 
eV for the first excited C4(^ng) state. Using Q and Mj conser- 
vation, we obtain that only the C(^Po) and C(^Pi, Mj - -1, +1) 
fine-structure substates are reactive; these three substates form 
vibrationally excited linear C4 in its fundamental state, while 
the other six substates lead to the formation of the ^Ilg excited 
state (through a high barrier). Any others approach on the triplet 
surface, as the perpendicular one leading to the rhombic C4, 
shows high barriers at the CASSCF/cc-pVDZ level. As the C4 
is formed with a great excess of energy, the potential minima 
of several electronic states, not just the ground state, lie below 
the dissociation limit (the potential energy of ground-state re- 
actants.). Stabilisation of C4 adduct can then occur via purely 
vibrational transitions from the ground-state adduct, but also via 
electronic transitions fro m excited state adducts if they can be 
formed without barriers dHerbst & Bates| [T988l). Here, we must 
therefore consider three mechanisms for adduct stabilisation: 

- vibrational radiative transitions from quasi-continuum vibra- 
tional levels in the ground state (X ^Eg) adduct to lower 
bound vibrational levels, 

- "inverse internal conversion" from the X ^E^ adduct to an ex- 
cited electronic state followed by electronic photon emission 
to bound vibrational levels of the ground state (case 1), 

- direct radiative emission from quasi-continuum vibrational 
levels in the ground X ^Eg state to lower vibrationally bound 
levels in excited states (case 2). 



State 


To (eV) th/exp 


R(Ci-C2) 


R(C2-C3) 


R(C3-C4) 


ei 


02 




0.89/0.93 


1.27 


1.33 


1.27 


180 


180 




0.63/0.82 


1.28 


1.34 


1.28 


180 


180 







1.32 


1.30 


1.32 


180 


180 



Distances R in A, angles in degrees and total energies are with respect 
to the X^Sg state. The results presented in the table are from this w ork 
except t he experimental data, which are from lLinnart et alj ( l200Clh for 
'Hg and lXu et alj ( 1 19971) for ^n,. 




d(C-C3) (1) 

Fig. 1. Co-linear VDZ-CASSCF one-dimensional cuts of the 
lowest triplet electronic states of linear C4 molecules along the 
C3-C stretch. The C3 bonds are optimized. The labeling of the 
C4 electronic states is based on D<x,h symmetry. 



Other possibilities arise from the isomerisation of the triplet 
linear C4 state toward the higher-energy singlet rombric C4 state. 
However this isomerisation has to overcome a high barrier as 
well as intersystem crossing. We can reasonably neglect those 
possibilities. 

Of the three itemized mechanisms, the case 2 situation can be 
ignored. Both experiments and detailed calculations on our part 
show that the potential surfaces of all relevant low-lying excited 
states of linear C4 lie almost always above the ground state at all 
relevant geometries so that emission transition from the ground 
state are impossible or very inefficient. Therefore, the only pos- 
sibility for stabilization via electronic emission is case 1, with 
the dominant channel being potential curve crossing from X ^Eg 

to ^n„ at long C2-C2 distances (see lMasso et al.l2006l for exam- 
ple) followed by an allowed u— >g electronic transition. The u-g 
potential curve crossing is not allowed when C4 lies in the sym- 
metric Dcoh linear configuration, but is allowed by vibrational 
excitation, which reduces the symmetry. Since the ^n„ state has 
an electronically allowed transition not only with the fundamen- 
tal state X ^Eg but also with the first electronically excited ^11^ 
state, we have to include these two transitions. 

If we designate unstable adducts ("complexes") with an as- 
terisk, the overall mechanism for association can now be written 
as: 



C-HC3 



(C4('Eg))* 



C + C3 



(C4('Er 



))• 



* C4(^E-) + hv 



(5) 



(6) 
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T(K) 

Fig. 2. RRKM lifetime (in s) of the C4 adduct as a function of 
the internal energy (in K) above the C + C3 dissociation limit. 



C4('n,)* — > 



C4{%,' Hg) + hv, 



(7) 



(8) 



where ass and diss designate association and dissociation, rad,v 
designates vibrational emission, cc designates curve crossing, 
rad,el designates electronic emission, and k designates a rate co- 
efficient. Assuming that curve crossing is more rapid than the 
assorted emission processes, we obtain the following expression 
for the rate coefficient of the overall radiative association (Icra): 

^RA = ^a,ss[^rad,v + ^lad.el-'^cc ] /(^rad.v + ^diss) , (9) 

where the equilibrium coefficient K^c is equal to fccc/^-cc and 
pertains to the ^II,, only . 

We first determine the microcanonical rate coefficient for 
dissociation, fediss, from the ground electronic state of linear C4 
as a function of internal energy E above the C + C3 dissocia- 
tion energy. We employ the so-called RRKM (Ri ce-Ramsperger, 
Kasse l, & Marcus) unimolecular approach ( Holbrook et al J 
11996 ^ using the FALLOFF system ( lForstlll993h . A so-called 
loose transition state is used with t he canonical flexible transi- 
tion state theory (CFTS) approach jHolbrooket al.|[T996 ). The 
needed structural and thermochemical parameters of the stud- 
ied system derive from quantum calculations. Thermal averag- 
ing of the microcanonical rate coefficient for dissociation of the 
C4 adduct obtained from the RRKM calculation allows us to cal- 
culate the thermal dissociation rate coefficient. The RRKM ther- 



mally averaged lifetime of the C4 adduct. 



1 



kji.s(T) 



is presented 



in Fig. 12] as a function of T. From the equilibrium constant gov- 
erning formation and dissociation of the adduct, the rate coeffi- 
cient for formation, which diff'ers from fcass in that fine-structure 
correlations are not yet taken into account, can be obtained. 

The calculated capture rate constant obtained by this proce- 
dure at 298 K is equal to 2.0 x 10"'" cm^^ s"' with very little tem- 
perature dependence in the 10-300 K range, in good agreement 
with theoretical treatments for reac tions between atom ic carbon 
and assorted alkenes and a lkynes ( Clary et al. 19941) and with 
the experimental results of Chast aing et al. for such reac- 

tions. This thermal rate coefficient, however, does not take into 
account the fraction of fine-structure states leading to the attrac- 
tive ground state C4 surface. With this factor and the assumption 



Table 5. Vibrational properties of the X ^E^, state of the C4 
molecule. 





V (cm ') 


IR int 

(km.mole"') 


Ai.o (s ') 


L ^Auo (s ) 


vl 


172 


27 


0.01 


0.34 


v2 


172 


27 


0.01 


0.34 


v3 


352 











v4 


352 











v5 


908 











v6 


1548 


312 


89 


342 


v7 


2032 












that the fine-structure states are in thermal equilibrium, we ob 
tain a refined association rate coefficient of f{T) x 2 x 10 
s ' in the range 10-1000 K, where 



10 cm^ 



f(T) = (1 + 2e-'^-4/T)/(i + 3 X e-'^-4/T + 5 x e-«-4/T). (lo) 

The method of estimating k,ad,v equate s this rate coefficient 
with the overall vibrational emission rate ( lHerbstlll982l: ISmithI 



A(£vib) 



EE 



(0 



(11) 



where A^^^^„._i, the Einstein A coefficient for mode ; sponta- 
neously emitting from level n,- to n,-l, is equal to ni x Aj'^y ), 

Pn.*(Evib) is the probability of finding mode / in level n,, and the 
sums are over all modes / and all energetically accessible quanta 
n,. The vibrational frequencies and Einstein coefficients, calcu- 
lated at the B3LYP/vtz level, are listed in Table |4] leading to 

^rad.v — 349 S 

To estimate the rate coefficient fci-ad.ei for stabilization by 
emission of electronic radiation, we use the results of ab initio 
calculations that the electronic potential curves of the C4(^n„), 
€4(^2^ ) and €4(^11^) states have very similar shapes. Therefore, 
we can approximately ignore nondiagonal Franck-Condon fac- 
tors. With the additional assumption that the diagonal factors are 
unity, we can write that < v'|v" >= Sy'y, so that the spontaneous 
electronic emission rate coefficients are given by the expression 
(esu-cgs units): 



rad.el 



(64n*e^al/3h)v^\^''\^, 



(12) 



where e is the electronic charge, v is the transition frequency, ao 
is the Bohr radius, and \f/'\ is the electronic transition moment. 
From the ^Hj, state, there are allowed transitions toward the 
and the ^11^ state ( lLarssonlll983l) . We obtain that ferad.ei « 6 x 10"* 
s- 1 . Finally, we obtain the equilibrium coefficient K^c by assum- 
ing that the relative populations of the ^II,,, ^Tlg, and states 
are equal to their relative densities of vibrational states at an in- 
ternal C4 energy of 5.21 eV, the experimental exothermicity of 
reaction (|4]l. This assumption leads to an equilibrium coefficient 
for the ^O,, state only of 0. 15, and an overall radiative stabi- 
lization rate A;rad,v + ^i-ad.el 

K,, of 9350 s-i. 
Once the partial rate coefficients are all calculated, we can 
obtain the total rate coefficient for radiative association, A:ra as a 
function of temperature. The results are shown in Figure [3l the 
value at 15 K is 3 X 10"'^ cm-' s"'. The value at 15 K if only 
vibrational stabilization is allowed is a much lower 1 x lO"'-' 
cm"* s"', which serves as a reasonable lower limit. The value 
obtained with the < v'|v" >= Sy'y assumption is certainly an 
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T(K) 

Fig. 3. Radiative association rate coefficient for C4 production as 
a function of temperature. 



extreme upper limit, because the diagonal Franck-Condon over- 
laps are certainly lower than unity, the electronic transition mo- 
ment, assumed to be independent of ener gy, may well dec rease 
with increasing energy for certain states ( Bove et aLll2002l) . and 
the internal conversion from the ground X ^2" of C4 to excited 
electronic states may not reach equilibrium. Estimation of the 
size these factors have on the association rate coefficient is not 
facile, but a reasonable guess is that they contribute to produce a 
factor of ~ 0.25, leading to a reduced value at 15 K of 8 x 10"'"' 
cm-' s"'. We adopt as an error limit an uncertainty of + 75%, 
allthough this error does not quite cover the extreme maximum 
and minimum values discussed above. From the RRKM cal- 
culations, we obtain a temperature dependence for A:ra of ap- 
proximately r~ ', in excellent agreement with the approach of 
iBates & HerbstI ( [T988i) for a system with two rotational degrees 
of freedom. We obtain finally a temperature-dependent associa- 
tion rate coeflicient;tRA(r) of (4.0 + 3.0)x 10"'4x(r/300)-i cm^ 
s ' in the range 10-300 K. Compared with t he radiative associ- 
ation rate coefficient of 1 x 10 '" cm^ s ' in lSmith et"!!] (|2004 
and OSU.2003, our new value is lower by a factor of ^ 100 at 10 
K. This discrepancy between old and new values is larger than 
the formerly assumed uncertainty of an order of magnitude. 



3.2. C + C5 

For the reaction between C in its ground state and linear C5 in 
its ground 'S^ state, the formation of the linear Ce adduct in its 
ground state is strongly exothermic and spin-allowed. This 
case is then very similar to the C -H C3 one, and to determine if 
the electronic curves are attractive along the entrance channel re- 
action pathway, we performed similar ab initio calculations. Our 
results are also similar, with no baiTier for the C(,i^l,~) electronic 
state but high barriers for the Cei^U) electronic states. Unlike the 
case of the C4 adduct, however, the Cg adduct can fragment by 
an exothermic route; in particular 



4 6 
log(time(yr)) 



Fig. 4. Average fractional difference (AFD, see text) as a func- 
tion of time. The dashed lines are for case (1), the dashed-dotted 
lines for case (2), and the solid lines for case (3). Black lines are 
for case (a) whereas grey curves are for case (b). 



to allow intersyste m crossing, as in the C + C2H2 reaction 
(iMebel et al.ll2007|) . and the baiTier on the singlet exit channel 
is certainly much lower than the exothermicity (' Nelson et al.l 
1982). The ffrst excited triplet state of C3 lies about 2 eV above 
the ground singlet state so is probably not accessible. Based 
on experiments, the assumption that the rate coefficient is in- 
dependent of the number of hydrogen atoms on the C5 reac- 
tant, and using the fine structure factor in eq. ( fTOl i , we estimate 
the rate coefficient to form the adduct and subsequent products 
to be 1.2 + (0.8) X 10"'° cm^ s"' for 10-300 K. This value is 
some what less than that obtained from the capture methodol- 
ogy of IClarv et al.l (|1994|) . Compared with the value in osu.2003 



-,-10 , 



,3 



C -H C5 is approximately the same but the products are C3 -1- C3 
instead of C(,. 

The rate coefficients for other C -H Cn reactions are discussed 
in the Appendix. These reactions are not among those the model 
results are very sensitive to. 



4. Impact on predicted abundances 

4. 1 . General impact 

To quantify the overall impact of change in the rate coeffi- 
cients/products computed in the previous section, we sum R'jit), 
the fractional difference for species j and reaction / as defined 
in equation ([TJ, over the species j and divide by the number of 
species nj to obtain the average fractional difference AFD(t) for 
reaction /: 



AFD(t) = i V \R'j(t)\. 
rij ^ ■' 



(14) 



Here, as in equation ([T]), we define the reference fractional abun- 
dance Xj-^(t) for each species as that obtained with the old rate 
coefficient/products, but define the new values of X'j{t) as those 
obtained with the modified values for reactions / obtained here. 
An unimportant reaction at a specific time would give an AFD 
close to zero at that time. It should be noted that although we 
C3('2-i-g)-i-C3('2-i-^e using a method similar to section 2.2, the goal here is com- 
(13) pletely different. In the sensitivity analysis section (2.2), we mul- 
with an overall reaction enthalpy of AH = -1.36 eV. Although tiply all reactions, one after another by a factor of 2 in order to 
this reaction is formally spin-forbidden to the products in their compare the abundance sensitivity to individual reactions. Here, 
singlet ground states, the adduct is likely to last long enough we change the rate coefficients for which new calculations have 



C('P)+C5('S+,) C6('I-,) C6('2) 
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Fig. 5. Abundances with respect to total hydrogen of several im- 
portant species. Black dash-dotted lines (0) are calculated abun- 
dances when the standard network osu. 2003 is used; grey dashed 
lines (1) when only C + C3 is modified; grey dash-dotted lines 
(2) when only C + C5 is modified, and black solid lines (3) when 
both C + C3 and C -t- C5 are modified. 



been done, for instance the C-1-C3 reaction is changed by two 
orders of magnitude, and we look how much the new rate co- 
efficients have changed the model predictions. If we redo the 
sensitivity analysis as in section 2.2 with the updated network, 
we may find different "important" reactions since the network 
has been changed. 

For the dense cloud model discussed in Section |2l we con- 
sider three different cases for reactions /: (1), we change only 
the rate coefficient of the C+C3 reaction, (2), we change only 
the products and the rate coefficient of the C+C5 reaction and 
(3), we perform the changes for both reactions. For each of these 
cases, we compute the AFD parameter by (a) summing over all 
species (452), and by (b) summing only over those species with 
a fractional abundance greater than 10"'" at 10^ yr For (b), the 
sum is made over 79, 80 and 77 species for cases (1), (2) and (3) 
respectively. Fig. |4] shows the results for the six different cases. 
It is clear from the results that the most important effect is ob- 
tained by the combined effect of changing both C+C3 and C+C5 
reactions. For case (3), the sum over all species yields an AFD 
of 15 at slightly more than 10^ yr, while the sum over the most 
abundant species only yields an AFD of 13 at 10^ yr If only 
one of the two reactions is changed, C-1-C3 has the greater effect 
whereas C-1-C5 affects the species at early times. 

4.2. Impact on specific species 

In Fig. |5] we show how the abundances of selected species are 
affected by the modifications of the C-1-C3 and C-i-Cs reactions at 



intermediate times; at other times there is little to no effect. For 
carbon monoxide, the reduction of the C -1- C3 reaction rate coef- 
ficient results in the lowering of its abundance in the period 10^ 
- 10^ yr by at most one order of magnitude. Changing the prod- 
ucts of the C+C5 reaction on the contrary slightly increases the 
CO abundance. For neutral atomic carbon, the gradual decrease 
in abundance with increasing time occurs much later when the 
rate coefficient of C +C3 is reduced; at its maximum, there is 
a difference in abundance of about 3 orders of magnitude. The 
change in C H- C5 speeds up the destruction of C. In fact, the 
only species shown in Fig. |5] for which the modification of this 
reaction has a strong effect are C5 and CeH. 

As discussed in Section 12.3! with the use of the osu. 2003 
network, these two reactions are both involved in catalytic cy- 
cles to convert C and O into CO. By radiative association re- 
actions with atomic carbon, C3 and C5 form C4 and Ce, which 
then react with atomic oxygen to give CO and regenerate C3 and 
C5. In this work, we have shown that the C H- C3 reaction is less 
efficient than the value used in osu.2003 by two orders of magni- 
tude at 10 K and the more likely products of the C H- C5 reaction 
are two C3 molecules. Even after changing these reactions, this 
catalytic process remains the most efficient one to form CO at 
lO"* yr. Changing only the C + C5 reaction products makes this 
process more efficient by introducing more C3 into the cycle. 
For this reason, CO is slightly increased whereas C is decreased. 
Obviously C3 is more abundant because C H- C3 consumes less 
C3 and C -H C5 produces it. Contrary to C3, the C5 abundance de- 
creases if the C -H C5 — » Ce - O -H Ce —> C5 + CO cycle is broken 
because Cg is not locked in this loop anymore and is destroyed 
to give C3. The effect of these changes on C-rich molecules such 
as the ones shown in Fig. |5]is less evident to explain. From a 
general point of view, larger abundances of atomic oxygen and 
carbon between 10"* and 10^ yr destroy these species more ef- 
ficiently than producing them. The radical CeH for instance is 
destroyed by C to give C7 and C7 in turn reacts with C to pro- 
duce C3 and C5 (see Table lA.ll . Similarly C4H2 is destroyed by 
the reaction C4H2 -H C ^ C5H + H and C5H is destroyed by O, 
producing CO and C4H. Large cyanopolyynes, such as HC7N, 
are reduced in abundance. On the other hand, large O-bearing 
molecules, such as CH3OH, are increased in abundance or not 
changed. 

4.3. Comparison witli observations 

To study the impact of these new rate coefficients on the agree- 
ment betwee n models and observat ions, we used the method 
described in IWakelam et al.l (l2006ah . We modeled the chem- 
ical evolution of two different clouds: L134N and TMC-1 
(Cyanopolyyne Peak), using the Nahoon chemical model. The 
list of reactions is osu.2003 with the original C + Cn rate co- 
efficients and OSU.2003 with the new values of the C -1- C„ rate 
coefficients listed in Table lA.ll For the uncertainties in the rate 
coefficients of these reactions, we assumed a flat distribution (be- 
cause the mean values of the rate coefficients are not necessar- 
ily the favorite values) with the uncertainty ranges described in 
Table lA.ll For the other reactions, we used log-normal distri- 
butions (see Wakelam et al. 2004). The physical conditions are 
taken to be homogeneous in the observed regions: a gas temper- 
ature of 10 K, an H2 density of lO"' cm"^ and a visual extinction 
of 10. The elemental abun dances are the low-meta l elemental 
abundances (see table 1 in Wakelam & Herbstll2008l ). For TMC- 
1 , we decreased the oxygen abundance to have a C/O elemental 
ratio of 1.2 in order t o better reproduce the observations (see 
IWakelam et al.l l2006ai and references therein). The observed 
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Fig. 6. Comparison between model and observations in two dif- 
ferent clouds (black: L134N and grey: TMC-1). The solid line 
represents the comparison with the model using the original val- 
ues of the C + Cn rate coefficients whereas the dashed lines rep- 
resents the comparisons with the updated model. 



abundances are given in Tables 3 and 4 o f iGarrod et al] ( l2007h . In 
addition, we included the observa tional upper l imits on the SiO 
abundances in these two clouds bv lZiurvs et alj(ll989l) : the upper 
limits for this species are 2.6 x 10'^^ and 2.4 x 10"^^ forL134N 
and TMC-1 respectively. In total, we have observed abundances 
of 42 and 53 species for these two clouds to compare with our 
calculations. 

The results of the comparisons are shown in Fig |6] in three 
panels, for each of which a different criterion is utilized. These 
criteria are based on the Monte Carlo approach to determina- 
tion of the uncertainties in the computed abundances based on 
uncertainties in the rate coefficients, and the overlap of Gaussian 
distributions representing both the computed and observed abun- 
dances (Wakelam et al. 2006a) . In the top panel, we use the frac- 
tion of reproduced molecules. In the middle panel, we use the 
so-called "distance of disagreement" D in which, for species 
not reproduced by the first criterion, we sum over the remain- 
ing species the distances between logarithms of the observed 
and model fractional abundances. Finally, in the bottom panel, 
we use the minimum distance of disagreement (over all times), 
divide by the distance of disagreement, and multiply by the 
fraction of reproduced molecules. See Wakelam et al. (2006^ 
for details on these parameters. The time of best agreement is 
roughly the same for all three parameters, although the sharpest 
agreement occurs for the bottom panel. It can be seen in all pan- 
els that the general agreement with the observations is slightly 
better using the new values of the C -i- Cn rate coefficients. 
However, the more interesting effect is the shift in time of the 
best agreement. With the original model, the best agreement was 



obtained around 3 x 10"^ yr whereas now we calculate the best 
agreement to be between 10^ and 3 x 10^ yr. N ote that to com- 
pare our results with the results of I Wakelam e t al. (2006a), we 
have only updated the C -t- Cn rate coefficients whereas many 
other reactions have been updated in the latest version of the osu 
database (osu.3.2008). 

5. Conclusion 

This paper presents an approach to improve chemical models, 
which consists of first identifying important reactions, then 
obtaining improved rate coefficients for these reactions, and 
finally studying the impact of these new rate coefficients on the 
predictions of the model. In this example of the method, we 
found that under cold dense cloud conditions, the abundances 
of many species are sensitive to the rate coefficients of the 
C -I- C3 and C -i- C5 reactions. Using detailed ab initio and 
kinetic calculations, we have shown that the rate coefficient of 
the C -H C3 — > C4 previously used in the osu. 2003 database is 
incorrect by two to three orders of magnitude at 10 K and has 
an inverse temperature dependence. We also found that the rate 
coefficient for the C -H C5 reaction was roughly correct in the 
osu. 2003 database, but the products are most likely C3 + C3 
instead of Ce, as assumed earlier based on faulty thermodynamic 
information. Newly obtained rate coefficients for other reactions 
in the C -H Cn series are also reported, although the model 
shows little sensitivity to these. With the new rate coefficient 
for the association reaction between C and C3 and the slightly 
changed rate coefficient and new products for the C -1- C5 
reaction, we found that the abundances of many carbon chain 
species are strongly modified at the so-called early time of 
10^ yr. The main reason for this effect is that the C H- C3 ^ C4 
/ O + C4 CO + C3 and C + C5 ^ Ce / O + C6 ^ CO -H C5 
reactions comprise a very efficient catalytic cycle to transform C 
and O into CO. Slowing down the first cycle and removing the 
second one, both decreases the abundance of CO and enhances 
the abundances of atomic C and O at early times. These changes 
in turn cause an alteration in the early-time abundances of many 
other species. The agreement with observed abundances in the 
two dense clouds L134N and TMC-1 (CP) is not significantly 
modified by these new rate coefficients/products but the best 
agreement is obtained at later times using the new values. 
The newly calculated rate coefficients and/or products will 
be included in the latest network from the Ohio State group 
(http : //www . physics . ohio- state . edu/~ eric/research. html). 
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Appendix A: Other C + Cn reactions 

The calculated rate coefficients and products for C + Cn are 
listed in Table lA.ll for n = 2 - 8 along with estimated uncer- 
tainties and the older values used in the osu.2003 network. The 
applicable temperature range is generally 10-300 K, although 
for radiative association reactions, the dependence can extend 
to higher temperatures. For the C(^P) + CiC'Sg) reaction, the 
spin-allowed pathway leads to a C3(a ^'llu) adduct, which is not 
the ground electronic state of C3 but is highly thermodynami- 
cally accessible, lying 5.16 eV below the reactants. Preliminary 
ab intio calculations (CASCCF, vdz basis, full valence active 
space) give no barrier for the production of the a^Ilu state. To 
determine the rate coefficient for radiative association, we con- 
sidered only vibrational relaxation. Computing vibrational fre- 
quencies for C3 at the B3LYP/vtz level, we obtained a relaxation 
rate of A;i-ad,v = 220 s"'. We then estimated the association (cap- 
ture) rate coefficient to be 2x 10"'" cm-' s"' based on experiments 
of analogous systems (e.g. C + C2H2) and calculated the disso- 
ciation rate coefficient by the RRKM method. We obtained the 
radiative association rate coefficient, shown in Table I A. ll to be 
(3 + 2) X 10""" X (r/300)"' " cm^ s"', a rather low value due to 
the rapid dissociation of the small C3 adduct. 

For the C(^P) + CaC^), C6(^2g), C^C^) reactions, the 
formation of the Cn+i adducts is strongly exothermic and spin- 
allowed. The possible potential energy surfaces in the en- 
trance channels for all of these reactions are and '•^•^H. 
Preliminary ab intio calculations (CASCCF, vdz basis) show that 
the three ''^'^1;+ states have no barrier in the entrance channel, 
but that the three '•^•^11 states have high barriers along the re- 
action pathway. Exothermic fragment channels exist for all of 
these adducts. The rate coefficients for these reactions are once 
again estimated from experimental values of analogous systems, 
which in general are smaller than the capture values, derived 
from the formalism of Clary et al. ( 1994), and further reduced 
by a factor that allows for the multiplicity of the potential en- 
ergy surfaces that correlate with the reagents and also by a factor 



Table A.l. Rate coefficients for the C + Cn reactions. 



Reaction 




osu.2003 rate coefficient 


New rate coefficient 














a 


P 


C + C, 


-^C3 


+ hy 


- 




- 


(3 ± 2) X 10-'" 


-1 


C + C3 


-^C4 


+ hy 


1.0 X 10" 


10 





(4 ± 3) X 10"'-* 


-1 


C + C4 


-^C2 


+ C3 


1.0 X 10" 


13 





(0.9 ± 0.6) X 10-'" 





C -1- Cs 


— > Cb 


+ hy 


1.0 X 10" 


10 









C + C5 




-1- f 

+ L3 








(1.2 + 0.8) X 10-'" 





C + Cs 


^C2 


+ C5 


1.0 X 10- 


13 





(0.5 ± 0.3) X 10-'" 





C + Cs 


^C3 


+ C4 








(0.5+ 0.3) X 10-'" 





C + Ci 


^Cg 


+ hy 


1.0 X 10^ 


10 









C + C7 


-^C3 


+ C5 








(1.4 ± 0.8) X 10-'" 





C + Ci 


-^C4 


+ C5 








(0.4 ± 0.2) X 10-'" 





C + Cg 


-^C3 


+ C6 








(0.4 ± 0.2) X 10-'" 





C + Cg 


— »C2 


+ Cv 


1.0 X 10- 


13 




(0.4 ± 0.2) X 10-'" 






' Rate coefficients in units of cm^ s ' given by A; = o- x (7/300)'' 



for the multiplicity of the exit channels if relevant. Correlation 
involving the fine-structure states of atomic C is ignored. The 
assigned rate coefficients and estimated uncertainties are shown 
in Table O 

The C(^P) + C7('Sg ) reaction leading to the formation of the 
Cg adduct in the ground electronic state is strongly exother- 
mic and spin-allowed. By comparison with C + Cg we assume 
that there is no barrier in the entrance channel for the formation 
of the ground C% electronic state and, although this reaction is 
formally spin-forbidden to the products in their singlet ground 
states, the adduct is once again likely to last long enough to al- 
low intersystem crossing. Once again, we use analogous expe ri- 
mental results rather than the formalism of lClarv et alj (1 19941) to 
obtain the capture rate, and reduce this via eq. ( flOl l. 



